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Abstract: The variation in soil organic carbon density (SOCD) has been widely documented at various spatial 
and temporal scales. However, an accurate method for examining the attribution of explanatory factors for 
change in SOCD is still lacking. This study aims to attribute and quantify the key climatic factors, anthropogenic 
activities, and soil properties associated with SOCD change in the native grasslands of Inner Mongolia, China, 
by comparing data between the 1960s and the 2010s. In 2007 and 2011, we resampled 142 soil profiles which 
were originally sampled during 1963—1964 in the native grasslands of Inner Mongolia. SOCD was determined 
in A horizon (eluvial horizon) of the soil. We selected the explanatory factors based on a random forest method, 
and explored the relationships between SOCD change and each of the explanatory factors using a linear mixed 
model. Our results indicated that the change in SOCD varied from the east to the west of Inner Mongolia, and 
SOCD was 18% lower in the 2010s than in the 1960s. The lower SOCD in the 2010s may primarily be 
attributed to the increasing in mean annual water surface evaporation, which explained approximately 10% and 
50% of the total variation and explainable variation in the change in SOCD, respectively. The sand content of 
the soil is also a significant explanatory factor for the decrease in SOCD, which explained about 4% and 21% of 
the total variation and explainable variation in the change in SOCD, respectively. Furthermore, the collection of 
quantitative information on grazing frequency and duration may also help to improve our understanding of the 
anthropogenic factors that govern the change in SOCD. 
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1 Introduction 


Soil is the largest carbon pool in the terrestrial biosphere, containing approximately 1550 Pg (1 
Pg=10!5 g) of soil organic carbon (SOC) and 950 Pg of soil inorganic carbon (Lal, 1999; Watson 
and Noble, 2002; Batjes, 2014). Any changes in SOC may affect the components involved in 
carbon cycling of terrestrial ecosystems by influencing the atmospheric CO2 concentration and 
soil-derived ecosystem services (Davidson and Janssens, 2006). Thus, SOC has significant 
impacts on the carbon balance in terrestrial ecosystems (Davidson and Janssens, 2006; Trumbore 
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and Czimcezik, 2008). 

Grasslands, covering about 30% of the global land surface and storing 28%—37% of carbon in 
terrestrial ecosystems, have been confirmed to explore the potentials for carbon sequestration (Lal, 
2004). As a gigantic carbon pool, grassland ecosystem plays a critical role in the global carbon 
cycling (Cao and Woodward, 1998). The grassland ecosystem is the largest terrestrial ecosystem 
in China (41% of the total land area), especially in Inner Mongolia Autonomous Region. Since the 
1960s, there have been enormous cultural changes in Inner Mongolia which have significantly 
impacted the grazing and land management in this region, leading to increasing pressures on the 
grassland ecosystem such as a rise in livestock numbers (increasing by 2.8x10° cattle and 3.1x10° 
sheep from the 1960s to the 2010s) and overgrazing. At the same time, regional climate has also 
changed in Inner Mongolia, with the annual mean temperature rising from 3.5°C to 5.1°C and 
mean annual precipitation decreasing from 320 to 272 mm (Hu et al., 2015). Such simultaneous 
changes (anthropogenic and climatic factors) can inevitably impact the SOC in the grassland 
ecosystem of Inner Mongolia. 

SOC storage is determined by the balance between the mass of carbon inputs and outputs, 
which are affected by many factors, such as climate change, land use change and soil properties 
(Burke et al., 1989; Bai et al., 2012; Butler-Lapointe, 2014; Zhao et al., 2015). For example, 
desertification could reduce SOC storage by decreasing the clay content in soil (Schjgnning et al., 
1999), as we know that clay may improve the bonds between soil organic matter (SOM) and soil 
particles, protecting soil aggregates and reducing SOC leaching. Precipitation events can promote 
soil respiration and increase net primary productivity in temperate grasslands, resulting in the 
accumulation of SOC storage (Chen et al., 2004; Yang et al., 2010). Increasing temperature has an 
acceleration effect on SOC decomposition, especially for stable SOC loss (Davidson and Janssens, 
2006; Lefèvre et al., 2014). Grazing can have profound impacts on soil properties and processes 
(Pifieiro et al., 2010). Since these above-mentioned factors may act simultaneously on SOC, it is 
difficult to quantify their relative contributions to the change in SOC at a regional scale and in a 
long-term scale. 

There are many researches related to the estimation of SOC storage and change (e.g., Li et al., 
2011; Dai et al., 2014; Conforti et al., 2016). However, the results may differ at regional scales or 
even at ecosystem scales (Goidts and Van Wesemael, 2007; Mestdagh et al., 2009; Fujisaki et al., 
2015; Conforti et al., 2016). For the grasslands of Inner Mongolia, Xie et al. (2007) estimated a 
loss of 0.56 Pg of SOC from the 1980s to the 2000s, while Dai et al. (2014) reported a slight 
increase of SOC storage in the topsoil (0.12 Pg C), with an increasing rate of 7 g C/(m?.a), over 
the same period. The native grasslands of Inner Mongolia, located on the Eurasian Continent, are 
an important area for climate change research along the Northeast China transect of the 
International Geosphere-Biosphere Program (IGBP), since this region is extremely sensitive to 
climate change (Dai et al., 2014). This region is also an important livestock production area and 
acts as a natural green protective barrier for North China (Zhang et al., 2011). In this study, we 
used two advanced methods (linear mixed model (LMM; also referred to as 
multilevel/hierarchical model) and random forest method) to analyze the attribution of factors 
impacting the change of soil organic carbon density (SOCD; a major indicator of SOC storage 
estimation) in the native grasslands of Inner Mongolia. LMM, which is gaining attention for use 
in ecological data analysis, is a generalized form of the standard linear model. In the generalized 
form, data are permitted to exhibit correlation and non-constant variability (Goidts et al., 2009). 
Furthermore, LMM estimates not only the means of the data (fixed effects) but also the variance 
and covariance (random effects). A random forest model was used to identify the subset of 
non-redundant explanatory factors introduced into the LMM. The random forest method is a 
popular tool for classification using high-dimensional data and can rank candidate predictors via 
inbuilt factor importance measures (Breiman, 2001; Janitza et al., 2016). 

Generally speaking, the overall objectives of this study were (1) to identify the long-term 
(approximately 50 years) effects of several key explanatory factors (including climate factors, 
anthropogenic factors and soil properties) on the change in SOCD and to (2) quantify the 
attributions of these explanatory factors on the change of SOCD in the native grasslands of Inner 
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Mongolia. 


2 Materials and methods 


2.1 Study area 


This study was conducted in the native grasslands of Inner Mongolia Autonomous Region 
(38°-S1°N, 107°-125°E; Fig. 1), China. Inner Mongolia has a temperate continental monsoon 
climate, with annual mean temperature ranging from —5°C to 9°C, mean annual precipitation 
ranging from 150 to 500 mm, and mean annual evaporation ranging from 1200 to 2500 mm, from 
the northeast to the southwest (Dai et al., 2014). 

The native grasslands of Inner Mongolia consist of three ecotypes, i.e., meadow steppe, typical 
steppe and desert steppe. The meadow steppe ecotype, located in the east of Inner Mongolia, is 
mainly distributed in the transitional zone between mountain coniferous forest and steppe around 
the Da Hinggan Mountains. The soil associated with this ecotype is predominantly chernozem and 
the vegetation is characterized by perennial mesophyte and xerophyte plants, with the dominate 
species of Stipa baicalensis, Filifolium sibiricum and Leymus chinensis (Ma et al., 2008). The 
typical steppe ecotype is located at the center of Inner Mongolia, where the soil is dominated by 
chestnut and the vegetation is composed of typical perennial xerophyte plants, with the 
constructive species of Stipa grandis, Stipa kryovii, L. chinensis and Stipa bungeana. The desert 
steppe ecotype is distributed in the west of Inner Mongolia and exhibits brown calcic soil. The 
vegetation is dominated by dwarf plants with high drought resistance, such as Stipa klemenzii, 
Stipa glareosa and Stipa breviflora (Ma et al., 2008; Dai et al., 2014). 


2.2 Soil sampling and analysis 


Soil data were derived from two integrated surveys conducted in the native grasslands of Inner 
Mongolia. The first was an integrated soil survey conducted by the Chinese Academy of Sciences 
during 1963—1964 (Inner Mongolia and Ningxia Integrated Survey Team, Chinese Academy of 
Sciences, 1978). The soil sampling sites were determined on the basis of typical grassland and 
soil types, with the final selection made by visual inspection in the field. The soil types included 
aeolian soil, solonetz, black soil, brown pedocal, chernozem, castanozem, dark loessial soil, 
desert soil, loess soil, fluvo-aquic soil, and sierozem. A total of 163 soil sampling sites were 
identified in this soil survey. Information on each soil sampling site was collected, including 
geographical location, grassland type, land use type, soil genetic horizons, SOM content, and 
particle size distribution. 

To identify the change in SOCD at temporal scales, we revisited the locations of the previous 
soil survey during the growing seasons of 2007 and 2011. It should be noted that we carefully 
compared each soil sampling site between the 1960s (the first soil survey in 1963—1964) and the 
2010s (the second soil survey in 2007 and 2011) in terms of geographical location (longitude and 
latitude), topographical characteristics (elevation and slope), land use type, vegetation information 
and soil profile descriptions (e.g., soil type, horizon name and thickness). Thus, we identified 142 
sampling sites that were matched with the original sampling sites regarding the profile descriptions, 
vegetation types, and land use types. Of these 142 soil sampling sites, only 122 sites were located 
in native grasslands in 1963; and, of these 122 soil sampling sites, 10 sites were subjected to tillage 
before the second soil survey in 2007. Finally, we selected 112 soil sampling sites that always 
remained within native grasslands during the period 1963—2011 in this study (Fig. 1). 

During the soil surveys, soil pits were excavated with a spade to the depth of 1 m at each 
sampling site, and the coordinates and descriptions of soil genetic horizons of the A horizon 
(eluvial horizon) were recorded. Approximately 4—6 soil layers were set in each soil sampling site 
on the basis of soil genetic horizons. Three soil samples were taken in each soil layer and then 
mixed for testing gravel content, soil texture and SOC content. For details, soil samples were 
air-dried, microscopic roots were removed by sieving, and fine roots were picked out by 
employing the electrostatic method. The gravel content of the soil was determined as the 
proportions of the soil skeleton with diameters in the ranges of 1-3, 3-5, and >5 mm. The gravel 
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fragment content >2 mm was calculated from the proportion of the soil skeleton and the 
corresponding diameter range using an exponential function. Clay, sand and silt contents from the 
A horizon collected in the 2010s were obtained with a particle size analyzer (LS 13-320, 
Beckman Coulter, USA). The SOC content in both reference periods (1.e., 1960s and 2010s) was 
determined by the K2Cr207-H2SO. wet oxidation method of Walkey-Black (Walkley and Black, 
1934). For soil bulk density analysis, soil samples were collected using metal cylinders (100 cm? 
in volume) at each soil sampling site, and then the samples were oven-dried at 105°C and 
weighed following the methods of Shi et al. (2012). 
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Fig. 1 Distribution of native grasslands and locations of soil sampling sites in Inner Mongolia Autonomous 
Region. It should be noted that soil sampling sites were classified into three geographical groups, i.e., 
northeastern (NE) Group, middle (M) Group and southwestern (SW) Group. 


2.3 Calculation of soil organic carbon density (SOCD) 


We calculated the SOCD based on the relationship between SOC content, soil bulk density and 
soil layer thickness (Eq. 1). 

SOCD =>". [(1-Rv)XBD;xSOC;xT;]/100, (1) 
where SOCD; is the soil organic carbon density of the i'* soil genetic horizon (kg C/m?); n is the 
number of soil genetic horizons in the soil profile; Ry; is the volume percentage of the fraction >2 
mm (%); BD; is the soil bulk density of the i® soil genetic horizon (g/cm); SOC; is the soil 
organic carbon content of the i soil genetic horizon (g/kg); and T; is the thickness of the i soil 
genetic horizon (cm). Rv: can be calculated from the weight percentage of gravel fragment (Rmi, %) 
using Equation 2: 

Ryi=RmiXBD/BD,yz, (2) 
where BD,+ is the specific gravity of the gravel fraction (g/cm). According to the method of Xie 
and Li (2012), we defined BD,; as 2.1 g/cm? in this study based on the average bulky density of 
gravel. 


2.4 Explanatory factors and data collection 


We firstly collected 17 explanatory factors based on the sampling site selection criteria, and then 
excluded 7 explanatory factors according to the results of a multi-collinearity test. Thus, a total of 
10 explanatory factors were selected in this study (Table 1). The original climate data were 
obtained from the National Meteorological Data Center of China (http://cdc.cma.gov.cn/home.do). 
We interpolated the climate data to 1 kmx1 km grids using ordinary kriging method. Considering 
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the lag effect of climate factors on the change in SOCD, we obtained a long-term climatic 
information in the form of 10-year averages for two periods (1954—1963 and 2000-2009). The 
changes in climate factors during the two periods were then calculated via subtraction. It should 
be noted that the climate factors were 10-year averages in our study because a previous research 
revealed that long-time periods (such as 20 or 30 years) had a negligible effect on SOC 
accumulation and modeling outcomes (Allen et al., 2013). 


Table 1 Explanatory factors used to explain the change in SOCD (soil organic carbon density) 


Classification Factor Description 


Change in mean annual precipitation (MAP) during the study period 


Climate AMAP (mm) (MAP 2o10¢-MAP 60s) 
Change in annual mean temperature (AMT) during the study period 
AAMT (°C 
CO (AMT 2010 AMT 5605) 
AATTS (°C) Change in mean annual cumulative temperature >5°C (ATTS) during the study period 
(ATT52010s-ATTS 1960s) 
Change in mean annual water surface evaporation (EVAP) during the study period 
AEVAP (mm) EVAP:010-EVAP 560) 
AMRH (%) Change in mean relative humidity (MRH) during the study period (MRH2010-MRH 1960s) 
Soil Clay (%) Clay content of the A horizon 
Sand (%) Sand content of the A horizon 
Silt (%) Silt content of the A horizon 
SOCD mean (kg C/m?) Mean values of natural logarithm-transformed SOCD 960s and SOCD 2910s 
Grazing ALD (DSE/hm’) Change in livestock density (LD) during the study period (LD 2910-LD 1960s) 


Note: DSE, dry sheep equivalent; A horizon, eluvial horizon. 


The estimation of changes in soil texture during the reference periods was not calculated in this 
study because these changes are small with a wide error range and have a negligible effect on the 
change in SOCD in deeper soil layers (Lie et al., 2012). The wide range of error in the changes 
of soil texture may occur because of the different determination methods applied in different 
periods. Therefore, the present study mainly focused on the contributions of different soil texture 
components to the change in SOCD. 

The number of livestock in each county of Inner Mongolia in the two 10-year periods was 
derived from the rural socioeconomic data collected during annual surveys (Inner Mongolia 
Autonomous Region Bureau of Statistics, 1949-2017). The livestock density (LD), calculated as 
the average quantity of livestock divided by the grassland area, was used to represent the mean 
grazing intensity in the two reference periods. The unit measurement of LD was DSE (dry sheep 
equivalent, i.e., the amount of feed required by a two-year-old, 50 kg Merino sheep to maintain its 
weight) per hectare. 


2.5 Statistical modeling 


The distribution of SOCD in both reference periods followed a lognormal distribution. Therefore, 
the SOCD data were transformed into a natural logarithm form (logSOCD), and the log change in 
SOCD between the two periods was then subtracted (denoted as AlogSOCD). Unless otherwise 
stated, all of the analyses were performed using R statistical software (R Core Team, 2017). The 
random forest method was used to identify the subset of non-redundant explanatory factors that 
had the greatest effect on the variability of the change in SOCD before the LMM was constructed. 
The random forest model can rapidly estimate the importance of each explanatory factor to the 
change in SOCD by reordering the value of each explanatory factor in turn. The explanatory 
factors with the greatest importance to the model can result in the greatest decrease in the 
accuracy of the model. The accuracy can be represented by the increase in the internal measure of 
the mean squared error (MSE) (Breiman, 2001). To obtain a robust ranking of the importance of 
the explanatory factors, we run the algorithm 500 times to account for the variability in 
importance rankings. 

Considering the homogeneous group effect and the autocorrelations between the sampling sites 
in the broad study area, we divided the sampling sites into three geographical groups 
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(northeastern (NE) Group, middle (M) Group and southwestern (SW) Group) according to the 
longitude and latitude of the sampling sites based on cluster analysis. The three geographical 
groups represented the northeast, middle and southwest of Inner Mongolia, respectively (Fig. 1). 

LMM could provide information on both random and fixed effects related to the means and 
trends across multiple hierarchical organizational levels (Serneels et al., 2007). Based on the 
importance ranking of the explanatory factors determined by the random forest method, we fitted 
an LMM equation for each geographical group (Eq. 3): 

AlogSOCD=Sand+AEVAP+AMAP+AAMT+ALD+ASOCD mean+(1|Groups). (3) 

The first six terms on the right side of Equation 3 are the fixed effects of the LMM and their 
definitions are shown in Table 1. It should be noted that although the AAMT (change in annual 
mean temperature (AMT)) had a low importance ranking according to the random forest analysis, 
we still used it in the model by considering the temperature sensitivity of SOM decomposition 
(Davidson and Janssens, 2006). Furthermore, to address the regression to the mean effect between 
the original SOCD and the change in SOCD, we used the mean SOCD value in the two reference 
periods (i.e., (SOCD 1960s+SOCD2010s)/2) in the LMM, instead of the initial value (Oldham, 1962). 
The term 1|Groups represents the random effect of the LMM (ie., the covariation of different 
geographical groups with change in SOCD). The LMM was fitted using a restricted maximum 
likelihood approach. 

In this study, we introduced the absolute contribution rate (ACR) and relative contribution rate 
(RCR) of the explanatory factors to quantify the interpretability of each explanatory factor in 
relation to the change in SOCD based on LMM regression. The ACR and RCR can be calculated 
by Equations 4 and 5, respectively. 


ACR, = = x 100%, (4) 
SS; 

RCR, = a x 100%, (5) 
SS p 


where ACR; is the contribution percentage of the i" explanatory factor to the total variation in the 
change in SOCD (%); RCR; is the contribution percentage of the i™ explanatory factor to the 
explainable variation in the change in SOCD (%); SSr is the sum of squares of the change in 
SOCD; SSz is the regression sum of squares of the LMM; and SS; is the partial sum of squares of 
the i" explanatory factor. We obtained the SS; from the decomposition of the regression sum of 
squares based on the analysis of variance (ANOVA). Considering the robustness of the LMM, this 
study mainly focused on the quantitative explanatory power of the main effect of each fixed effect 
factor on the change in SOCD. In other words, interaction effects were not included in this study. 


3 Results 


3.1 Change in SOCD 


The SOCD at the depth of 0-1 m in the native grasslands of Inner Mongolia decreased from 9.95 
(+7.26) kg C/m? in the 1960s to 8.12 (+7.36) kg C/m? in the 2010s (decreased by 18%). The 
change in SOCD exhibited a decreasing trend from the northeast to the southwest of Inner 
Mongolia, especially in the middle part (Fig. 2). The soil bulk density increased from 1.38 (+0.04) 
g/cm? in the 1960s to 1.46 (+0.05) g/cm? in the 2010s. The SW Group presented the largest 
increase of soil bulk density (0.14 g/cm; Table 2). The spatial trend of change in SOC content 
was similar to that in SOCD. The SOC content was 0.46% lower in the 2010s than in the 1960s. 
The differences in the changes of SOC content and bulk density indicated that the decrease in 
SOCD in the M Group was mainly related to the decrease of SOC content, while the decrease in 
SOCD in the SW Group was resulted from the differences in bulk density. 


3.2 Explanatory factor selection 


The rank of importance of the 10 explanatory factors returned by the random forest analysis is 
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Fig. 2 SOCD (soil organic carbon density) in the 1960s (a) and 2010s (b) as well as the change in SOCD from 
the 1960s to 2010s (c) for each geographical group in the native grasslands of Inner Mongolia. The whiskers of 
the box and whisker plots denote the maximum value (less than or equal to the 3" quartile+1.5IQR (interquartile 
range)) and the minimum value (greater than or equal to the 1* quartile-1.5IQR), respectively. The circle points 
represent the values that greater than the 3" quartile+1.5IQR or less than the 1* quartile-1.5IQR. 


Table 2 Changes in SOC (soil organic carbon) content and soil bulk density at the depth of 0-1 m in the native 
grasslands of Inner Mongolia from the 1960s to the 2010s 


Change in SOC content Change in soil bulk density 


Group 


Change (%) Change rate (%/100 a) Change (g/cm?) Change rate (g/(cm?-100 a)) 
NE 0.09 0.20 0.02 0.04 
M -1.37 —2.80 0.07 0.14 
SW —0.09 —0.20 0.14 0.28 
Mean —0.46 —0.90 0.08 0.15 


Note: Soil sampling sites were classified into three geographical groups, i.e., northeastern (NE) Group, middle (M) Group and 
southwestern (SW) Group, representing the northeast, middle and southwest of Inner Mongolia, respectively. 


shown in Figure 3. Six explanatory factors can be considered to influence AlogSOCD: AEVAP 
(change in mean annual water surface evaporation (EVAP)), AMAP (change in mean annual 
precipitation (MAP)), sand content, SOCDmean (mean value of natural logarithm-transformed 
SOCD 960s and SOCD2010;), ALD (change in livestock density (LD)) and AAMT. Again, although 
the AAMT had a low importance ranking according to the random forest analysis, we used it in 
the model by considering the temperature sensitivity of SOM decomposition. 

The means and ranges of the six selected explanatory factors in each geographical group are 
shown in Table 3. Values of factors that correlated with temperature (AEVAP and AAMT) and 
grazing (ALD) increased in the 2010s compared with those in the 1960s, whereas AMAP 
decreased markedly over the past 50 years. SOCD mean and ALD exhibited the same decreasing 
trend as the change in SOCD from the northeast to the southwest of Inner Mongolia. Conversely, 
AMAP displayed an opposite spatial trend, i.e., showing an increasing trend from the northeast to 
the southwest. 


3.3 Contributions of explanatory factors to the change in SOCD 


Coefficients and variances for the random and fixed effects of the LMM are shown in Table 4. 
The variance of random effects indicated that 10% of the total variability in the change in SOCD 
existed between geographical groups. Therefore, we concluded that it is appropriate to use an 
LMM to analyze the change in SOCD. The change in SOCD was negatively correlated with 
AAMT and AEVAP and positively correlated with AMAP. The change in SOCD was significantly 
influenced by AEVAP and soil texture. 

Based on the LMM, we obtained the contributions of the explanatory factors to the change in 
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SOCD using Equations 4 and 5 (Table 5). Climate factors (AEVAP, AMAP, and AAMT) explained 
9.84% and 51.57% of the total variation and explainable variation in the change in SOCD, 
respectively. Among the climate factors, AEVAP contributed the most to the change in SOCD 
whereas AAMT contributed the least. Additionally, soil texture accounted for 3.94% and 20.63% 
of the total variation and explainable variation in the change in SOCD, respectively. The 
contribution of the variance explained by the grazing intensity was much lower than those 
explained by the other explanatory factors. 


AEVAP 4 t----{ | }----4 
AMAP 4 t----[ | }---4 


Explanatory factor 


AMRH 4 H----- | | }----4 


T T T T 


T 
2 4 6 8 10 12 14 
Increase in MSE (%) 


Fig. 3 Relative importance of the explanatory factors to the variation in the change in SOCD. The explanatory 
factors with the higher importance to the model can result in the higher decrease in the accuracy of the model 
(expressed as MSE (mean squared error)). The whiskers of the box-and-whisker plots denote data extremes and 
the importance of each factor decreases down the plot. AEVAP, change in mean annual water surface evaporation; 
AMAP, change in mean annual precipitation; Sand, sand content of the A horizon (eluvial horizon); Clay, clay 
content of the A horizon; SOCDmean, mean values of natural logarithm-transformed SOCD 1960s and SOCD2010s; 
ALD, change in livestock density; AMRH, change in mean relative humidity; Silt, silt content of the A horizon; 
AATTS5, change in mean annual cumulative temperature >5°C; AAMT, change in annual mean temperature. 


Table 3 Means and ranges of the selected explanatory factors in each geographical group in the native 
grasslands of Inner Mongolia 


Group Number Sand AMAP AEVAP AAMT ALD SOCD mean 
of sites (%) (mm) (mm) (°C) (DSE/hm?) (kg/m?) 
NE 40 53.63 -52 147 1.53 0.81 14.21 
(23.23-83.23) (125-37) (134-406) (0.06-2.72) (0.15-2.62)  (2.21-35.55) 
M 34 56.81 —60 138 1.68 0.33 8.49 
(11.97-83.23) (99-23) (9-278) (0.96-2.33) (0.12-1.17)  (1.05-24.36) 
sw 38 62.12 —10 142 1.93 0.19 5.10 


(12.22-85.97) (86-148) (326-531) (0.45-5.27) (0.22-0.80)  (0.97-15.46) 


Note: Means are listed in the front of the brackets and ranges are listed in the brackets. 
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Table 4 Random and fixed effects of the linear mixed model for the natural log-transformed change in SOCD 
(AlogSOCD) 


Random effects 


Groups Name Variance Standard deviation 
Geographical group Intercept 0.022 0.148 
Residual 0.189 0.435 
Fixed effects 
Parameter Estimate Standard error t value P value 
Intercept 0.284 0.252 —0.026 0.980 
SOCD mean 0.046 0.071 0.641 0.523 
Sand —0.004 0.002 -1.715 0.089 
AMAP 0.001 0.001 —0.212 0.833 
AAMT —0.030 0.067 0.166 0.869 
AEVAP —0.001 0.000 -3.502 0.001 
ALD 0.137 0.115 0.498 0.619 


Table 5 Contributions of the explanatory factors to the change in SOCD 


Explanatory factor SS; SSr SSr ACR (%) RCR (%) 
SOCD mean 0.12 0.48 2.53 
Sand 0.98 3.94 20.63 
AMAP 0.07 0.28 1.47 
4.75 24.90 

AAMT 0.02 0.08 0.42 
AEVAP 2.36 9.48 49.68 
ALSD 0.05 0.20 1.05 


Note: SS;, the partial sum of squares of the i explanatory factor; SSr, the regression sum of squares of the LMM; SS, the sum of 
squares of the change in SOCD; ACR, absolute contribution rate; RCR, relative contribution rate. 


4 Discussion 


Our findings revealed a SOCD loss of 1.83 (+4.57) kg C/m? in the native grasslands of Inner 
Mongolia during 1963-2011. This result was comparable to some previous observations in Inner 
Mongolia. For example, Wang et al. (2003) estimated a decreasing trend of SOCD from the 1960s 
to the 1980s in the grasslands of Inner Mongolia (from 2.15—5.00 to 0.00-2.15 kg C/m?) based on 
data from the first and second national soil surveys of China. Xie et al. (2007) detected a decrease 
of SOCD (decreased by 0.71 kg C/m?) in the grasslands of Inner Mongolia from the 1980s to the 
2000s using data from the Second National Soil Survey of China. These results indicated that the 
SOCD in the grasslands of Inner Mongolia decreased rapidly from the 1960s to the 1980s and 
decreased slowly from the 1980s to the 2000s. 

In this study, we detected a significant negative correlation between the change in SOCD and 
the original SOCD (r= —0.21, P<0.03), indicating that the carbon-rich soils tend to lose carbon. 
On the contrary, the carbon-poor soils generally tend to gain carbon. This phenomenon was 
identified through regression to the mean effect, which is a methodological bias that commonly 
occurs between the change and original values (Pearson, 1896; Chamberlain et al., 2010). Oldham 
(1962) suggested that the means of two measurements might be used, instead of the original 
values, to account for the fact that the samples were characterized by non-equivalent time 
intervals between two surveys (Lark et al., 2006). In this study, we introduced the natural 
logarithm-transformed SOCD mean in the LMM to identify the contribution of the original SOCD 
to the change in SOCD. Our results indicated that the original SOCD had a little effect on the 
change in SOCD, which explained only 0.48% and 2.53% of the total variation and explainable 
variation, respectively (Table 5). Thus, the significant negative correlation between the original 
SOCD and the change in SOCD is mainly due to the regression of the mean effect. In this study, 
we introduced an LMM to minimize the bias generated by spatial heterogeneity. To obtain a 
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robust regression model, we mainly focused the LMM on the direct effects of the explanatory 
factors on the change in SOCD and did not consider the interactions or indirect effects of these 
factors, since insufficient data were available to support such a complex model. It should be noted 
that we used three geographical groups derived from the cluster analysis instead of grassland 
types as random effects in the LMM, because the classification of grassland types is a 
comprehensive result considering the climatic, geographical, biological, and ecological 
characteristics. Therefore, grassland types are not independent of the explanatory factors and 
cannot be the random effects in the LMM. The geographical groups are derived on the theory of 
spatial autocorrelation. The relationship between two sampling sites is correlated only with the 
distance between them, i.e., the closer the two sampling sites, the more similar the two sampling 
sites. Thus, the geographical groups are independent of the fixed factors. 

Numerous studies have documented that climate change significantly influences SOC content 
at different spatial and temporal scales (Cao and Woodward, 1998; Fantappié et al., 2011; 
Shangguan et al., 2012; Allen et al., 2013). For example, Fantappié et al. (2011) proposed that a 
moderate or large decrease in MAP (decrease of >2 mm) and a moderate-to-large increase in 
AMT (increase of >0.62°C) generally had a stronger effect on the decrease in SOC content in 
meadows than in forests or arable land in Italy from 1961 to 2008. Climate factors related to soil 
moisture, such as precipitation and evaporation, which incorporate approximated the first-order 
factors, may be strongly related to the change in SOCD at a regional scale (Wynn et al., 2006; 
Allen et al., 2013). This study adds further evidence, as the direct effects of all explanatory factors 
explained about 15% of the total variation in the change in SOCD. The change in evaporation was 
the dominant explanatory factor on the change in SOCD, explaining 9.48% and 49.68% of the 
total variation and explainable variation, respectively (Table 5). The availability of water is 
mainly dependent on the rates of water input (precipitation) and output (evaporation). The amount 
of soil moisture output increases as evaporation increases. Then, the amount of water available for 
plant growth is limited and the rate of biomass input to soil is reduced, consequently decreasing 
the carbon flow transformed to the soil carbon pool (Wynn et al., 2006). In this study, the 
observed climate change-dominated SOCD dynamics should be considered spatially, as the 
sampling sites were distributed across a wide spatial range in Inner Mongolia. The decreases of 
SOCD in NE, M and SW groups were mainly driven by the changes in EVAP, MAP and AMT, 
respectively. These findings indicated that the decrease in SOCD dominated by climate factors 
could vary at regional scales. 

As shown in Table 5, the contribution of the sand content to the change in SOCD ranked the 
second among the selected explanatory factors, explaining 3.94% and 20.63% of the total 
variation and explainable variation, respectively. According to our study, SOCD exhibited a 
decreasing trend in sandy soil and an increasing trend in clayey soil in the grasslands of Inner 
Mongolia. These findings verified that fine soil texture tends to protect organic matter while 
coarse soil texture tends to loss organic matter (Burke et al., 1989; Parton et al., 1993, 1995). 
Silver et al. (2000) demonstrated that the flow of carbon from slow-turnover SOM to passive 
SOM is inversely related to the sand content, and Wiesmeier et al. (2015) concluded that the 
capacity of soil carbon storage in semi-arid steppes depends on the silt and clay contents. 

Issues regarding the effect of grazing on the SOC content are debated and no consensus has yet 
been reached (Frank et al., 1995; Reeder et al., 2004; Goidts et al., 2009; He et al., 2011). 
Generally speaking, grazing influences the SOC in a complex manner (Zou et al., 2007; Piñeiro et 
al., 2010; Butler-Lapointe, 2014; Pringle et al., 2014). In this study, the grazing density was found 
to have a little positive effect on the decrease in SOCD, which is in accordance with the findings 
of some previous studies conducted in Inner Mongolia (e.g., Zhao et al., 2007; Han et al., 2008; 
Steffens et al., 2008; Wiesmeier et al., 2009; Liu et al., 2012). Wiesmeier et al. (2009) determined 
the SOC stocks at two continuously grazed sites and two ungrazed (fenced and excluded) sites in 
steppe ecosystems dominated by L. chinensis and S. grandis in Inner Mongolia, and found that 
SOC stocks were significantly lower while bulk density was significantly higher in the topsoil 
layer (0—4 cm) at both continuously grazed sites. Steffens et al. (2008) sampled five sites with 
different grazing intensities (ranging from ungrazed site to heavily grazed site) at a depth of 0—4 
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cm in a semi-arid steppe of Inner Mongolia, and found a significant lower SOC stock in the 
heavily grazed site while no significant differences in the remaining sites. In this study, ALD was 
observed to be positively correlated with the change in SOCD, which may be due to the large 
uncertainty associated with the determination of livestock grazing density or the strong 
correlations of the change in SOCD with climate factors and the original SOCD (Chan et al., 
2010). 


5 Conclusions 


Based on two integrated soil survey datasets in the native grasslands of Inner Mongolia, we found 
a significant decrease of SOCD from the 1960s to the 2010s. The contributions of the explanatory 
factors (climate factors, soil texture, and grazing intensity) to the change in SOCD were 
quantified using a random forest model and a linear mixed model. The climate factors dominated 
the change in SOCD, especially the mean annual water surface evaporation, which explained 
approximately 10% and 50% of the total variation and explainable variation, respectively. Thus, 
climate may be the most critical factor in driving the SOCD dynamics in the temperate grasslands. 
Soil texture was the secondary controlling factor on the change in SOCD, especially the sand 
content, which explained about 4% and 21% of the total variation and explainable variation, 
respectively. The influence of grazing intensity in the change in SOCD was difficult to determine 
and only a small contribution was detected. To sum up, the results of this study have a profound 
significance for ecosystem management and conservation in the native grasslands of Inner 
Mongolia. 
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